#ifndef UHM_MAT_OPS_HXX_
#define UHM_MAT_OPS_HXX_

namespace uhm {

  // Mat_ Ops
  // --------
  void random_matrix(Mat_Base_ &A);
  void random_spd_matrix(Mat_Base_ &A);
  
  void identity_matrix(Mat_Base_ &A);
  void norm1(Mat_Base_ &A, Double_ &val);
  void norm_inf(Mat_Base_ &A, Double_ &val);
  
  void copy(Mat_Base_ &A, Mat_Base_ &B);
  void copyt(Mat_Base_ &A, Mat_Base_ &B);
  void axpy(Int_ alpha, Mat_Base_ &A, Mat_Base_ &B);


  
  // Interface
  // ---------
  Int_ apply_diag_scal     ( Int_ side, Int_ conj, Int_ uplo,
                             Mat_Base_ &A );
  Int_ apply_diag_scal     ( Int_ side, Int_ conj, 
                             Mat_Base_ &A, Mat_Base_ &B );
  Int_ apply_diag_inv_scal ( Int_ side, Int_ conj, 
                             Mat_Base_ &A, Mat_Base_ &B );

  Int_ apply_pivots( Int_ side, Int_ trans, Mat_Base_ &piv, Mat_Base_ &A );

  Int_ chol_nopiv  ( Int_ uplo, Mat_Base_ &A );
  Int_ chol_piv    ( Int_ uplo, Mat_Base_ &A, Mat_Base_ &p );
  Int_ chol_incpiv ( Int_ uplo, Mat_Base_ &A, Mat_Base_ &p );

  Int_ ldl_nopiv   ( Int_ uplo, Mat_Base_ &A );
  Int_ ldl_piv     ( Int_ uplo, Mat_Base_ &A, Mat_Base_ &p );
  Int_ ldl_incpiv  ( Int_ uplo, Mat_Base_ &A, Mat_Base_ &p );

  Int_ lu_nopiv    ( Mat_Base_ &A );
  Int_ lu_incpiv   ( Mat_Base_ &A, Mat_Base_ &p );
  Int_ lu_piv      ( Mat_Base_ &A, Mat_Base_ &p );

  Int_ trmm        ( Int_ side, Int_ uplo, Int_ trans, Int_ diag,
                     Int_ alpha, Mat_Base_ &A, Mat_Base_ &B );
  Int_ trsm        ( Int_ side, Int_ uplo, Int_ trans,
                     Int_ diag, Int_ alpha, Mat_Base_ &A, Mat_Base_ &B );
  Int_ gemm        ( Int_ transa, Int_ transb, Int_ alpha,
                     Mat_Base_ &A, Mat_Base_ &B, Int_ beta, Mat_Base_ &C );
  Int_ gemm_l      ( Int_ transa, Int_ transb, Int_ alpha,
                     Mat_Base_ &A, Mat_Base_ &B, Int_ beta, Mat_Base_ &C );
  Int_ gemm_u      ( Int_ transa, Int_ transb, Int_ alpha,
                     Mat_Base_ &A, Mat_Base_ &B, Int_ beta, Mat_Base_ &C );
  Int_ syrk        ( Int_ uplo, Int_ trans, Int_ alpha,
                     Mat_Base_ &A, Int_ beta,  Mat_Base_ &C );
  Int_ herk        ( Int_ uplo, Int_ trans, Int_ alpha,
                     Mat_Base_ &A, Int_ beta,  Mat_Base_ &C );

  // UHM_ Ops
  // --------
  void random_matrix(UHM_Base_ &U);
  void random_spd_matrix(UHM_Base_ &U);

  void identity_matrix(UHM_Base_ &U);

  void schur_x2r(UHM_Base_ &U);

  void copy(UHM_Base_ &A, UHM_Base_ &B);
  void copy_A(UHM_Base_ &A, UHM_Base_ &B);
  void copy_B(UHM_Base_ &A, UHM_Base_ &B);

  void axpy(UHM_Base_ &A, UHM_Base_ &B);
  void axpy_A(UHM_Base_ &A, UHM_Base_ &B);
  void axpy_B(UHM_Base_ &A, UHM_Base_ &B);
}


#endif
